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Abstract 

The lowest-lying collective modes of a trapped Bose gas in an optical lattice are studied in 
the Bose-Hubbard model. An exact diagonalization of the Hamiltonian is performed in a one- 
dimensional five-particle system in order to find the lowest few eigenstates. Dipole and breathing 
character of the eigenstates is confirmed in the limit where the tunneling dominates the dynamics, 
but under Mott-like conditions the excitations do not correspond to oscillatory modes. 
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I. INTRODUCTION 



Exciting collective modes is a useful and popular tool for probing the many-body physics 
of trapped atomic gases. Following the first creation of a trapped condensate in 1995 1], 
modes in these systems have been subject to extensive theoretical and experimental study |2|. 
The fundamental zero-temperature theory was laid down about half a century ago [jj, and 
was readily adapted to the case of trapped condensates In three dimensions, oscillatory 
modes are naturally classified according to their multipolarity, and can be selectively excited 
by deforming the magnetic trap, or applying laser pulses that repel or attract the atoms in 
selected regions of space 0, f|. 

The picture is complicated considerably if one adds an optical lattice, consisting of one 
or several standing laser waves that act as a spatially periodic potential on the atoms. In 
the limit of a weak optical potential, the mode frequencies are simply given by those of the 
trapped Bose-Einstein condensed cloud in the absence of an optical lattice, but renormalized 
by the effective mass acquired by the bosons in the periodic potential p, 0]. Away from 
this limit, however, the presence of an optical lattice offers quite different physics, and a 
new phase appears, namely the Mott insulator, when the interactions are strong HQ. In 
addition, when an external trapping potential is present, there exist parameter regimes where 
spatially separated regions of Mott-insulating and superfluid behavior coexist The 
behavior of the trapped gas and the nature of its collective modes are expected to become 
quite different in these regimes compared to the quite well understood case of a trapped 
cloud with no optical lattice present 0|. 

In order to be able to address both the strongly and weakly interacting case and the 
crossover between these, we shall study the Bose-Hubbard model in the exactly solvable 
case of few bosons and in one dimension. The method shall be exact diagonalization in 
a truncated basis. This way we hope to gain qualitative knowledge of the spectrum that 
applies also to larger systems and higher dimensions. The paper is organized as follows. The 
Hamiltonian and the numerical method are explained in Sec. |TI] The nature of the ground 
state and the low-lying excitations in a shallow trap is investigated in Sec. II 111 and the case 
of a tight trap in Sec. IIVI Concluding remarks are given in Sec. El 
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II. BOSE-HUBBARD MODEL AND TRUNCATED BASIS 



The starting point is the Bose-Hubbard Hamiltonian 12 



H — —a\.a\,a r a r — — a\(a r+ i + a r _i) H — aja r . (1) 

The first term in this Hamiltonian describes the interactions which are effectively repul- 
sive if U > (which is the case in this paper), the second, so-called tunneling or hopping 
term is associated with the kinetic energy, and the last term describes the external trap- 
ping potential. The index r denotes the spatial position and takes on integer values. Such 
a one-dimensional Hubbard model describes a system with a tight trap in the directions 
perpendicular to the lattice so that the other degrees of freedom are frozen out, thus re- 
sembling a coupled chain of quantum dots. Higher dimensions will make the picture more 
complicated, but the main qualitative features observed in the present paper are expected 
to carry over to higher dimensions. 

The Hamiltonian contains three physical parameters. The tunneling strength J can be 
written 

m or z 

where 8r is the spacing between wells and m* is the effective mass acquired by the atoms 
due to the periodic potential 0. The interaction strength U is related to trap parameters 



through the relation 



U 



4irh 2 a 



J d 3 r|v|/ T B(r)| 4 , (3) 



m 

where a is the s-wave scattering length and ^xb is the ground-state wave function in one 
potential well in the tight-binding approximation. The effective trap frequency u is defined 
in terms of the bare particle mass m and trap frequency Q as 

uj = ^/mVL8r. (4) 

Let us at this point rescale the Hamiltonian and work in units of J; formally we set J = 1 and 
retain U and uj as the two parameters of the system. In addition to the parameters already 
discussed, the number of atoms iV or alternatively the chemical potential is a parameter of 
the system; we shall fix iV = 5 in this study. Furthermore, for a few-particle system the 
even/odd parity of the number of sites L may also play a decisive role; such effects vanish in 
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the limit of large systems. For definiteness only odd L will be considered, but some attention 
will be paid to parity effects where appropriate. 

In a number-conserving formalism, the natural basis is the set of real-space Fock states 
that are also eigenstates of the interaction and trap energies: 

| . . .n r _ 1 n P n r+1 ...) = ... (4-i) n - 1 (4) nr (oJ+i) nH - 1 • • • l°>- ( 5 ) 

The trapping potential implies a finite system size: it turns out that between one and 25 
sites is needed to accommodate a system of five particles for the trap parameters considered 
here. The size of the basis for a system of size L with N particles is (N+L — 1)!/[JV!(L — 1)!]; 
for iV = 5, L = 25 the number of states is 118755. Clearly, the computations can be made 
much more efficient if the basis is truncated so that the many improbable Fock states do not 
contribute: it is immediately obvious that states such as |iV000 ... 0), where all the particles 
are concentrated at one endpoint of the lattice, make only a very small contribution to the 
dynamics. 

There is, therefore, much to be gained if the basis is truncated. The following scheme 
turns out to be practical for both strong and weak coupling, although it was designed for 
dealing with Mott-insulator-like conditions where the interactions are strong. Start with a 
single Fock state labeled |1'), for instance the state with all the particles at the same site, 
|1') = | . . . OOiVOO . . .). Now enumerate all the states, |2'), . . . , \n'^) y that can be constructed 
from by one application of the tunneling term in the Hamiltonian (fl|). and construct the 
Hamiltonian matrix elements in the process. Operate again with the tunneling Hamiltonian 
on the states |2'), . . . , |n^) to form new states \n[ + 1), . . . , \n' 2 ), taking care not to double- 
count states; iterate this step p times so that a basis is formed that consists of n' states. 
Within this basis, the Hamiltonian is now diagonalized and the ground state is found. 
Among the n' p Fock states in the preliminary basis, choose the one that has the largest 
overlap with the ground-state eigenvector and label it Now discard all the other Fock 
states that were just constructed, and instead iterate the whole scheme again to construct 
a new basis |1"), . . . , \n'p); do the iteration a few (say, M = 3 or 4) times. The basis thus 
constructed, . . . , |?tX m )), will contain all the Fock states that have significant overlap 

with the ground state for the given physical parameters; in a sense, by constructing this 
basis tunneling effects to pth order have been incorporated. Convergence with respect to 
p and M is readily checked, so that the diagonalization can for all practical purposes be 
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considered exact. 

The diagonalization is performed with ARPACK, which uses an Arnoldi algorithm. 



III. MODES IN SHALLOW TRAPS 

The competition between the tunneling, interaction and trap energies gives rise to a 
rich phase diagram (cf. 0, Q]). Consider first the shallow trap. Figure Q displays the 
ground-state density distribution for values of U ranging from weak to strong interactions, 
with u = 0.3. The quantum fluctuations of the number of particles in the central well are 
displayed in Fig. |21 It is seen that we are in the fluctuation-dominated, superfluid regime. 
The very slight suppression of fluctuations in the centermost well is in fact a signal that we 
are in the vicinity of the Mott insulating regime; if the trapping strength is increased, this 
suppression becomes stronger, as shall be discussed in Sec. I1VI For the present trapping 
strength, u> = 0.3, the effect is barely noticeable even for U = 100. 

In this shallow trap, the competition between the tunneling and interaction is decisive in 
the interior of the system, while the trap still determines the size and the behavior close to 
the boundary. The density is spread out over many sites, and therefore the system resembles 
a trapped, Bose-condensed cloud in the absence of a lattice. It is in this "trapped BEC" 



regime that the experiments of Ref. [7] were conducted. While five particles are too few to 
be considered truly in the trapped BEC regime, the results reported here still give a hint of 
that limit, as is seen in Fig. ^ In a trapped condensate in the absence of an optical lattice, 
the density distribution is Gaussian for weak coupling, but flattens out and takes on the 



shape of an inverted parabola for stronger coupling |J. The density profile* shown in Fig. 
[T]can of course not be expected to exactly follow this behavior, because of the discreteness, 
but the dependence on coupling is similar. 

The mode frequencies, i. e. the excitation energies relative to the ground state, for the 
lowest four excited states are displayed in Fig. El When the coupling becomes weak, the mode 
frequencies approach integer multiples of the trap frequency with the expected degeneracies 
in a harmonic trap. Restoring dimensions, with the aid of Eqs. (I2I4D . we find in fact 0,0] 



/ TYl 

E n = qyj—Ml, (6) 
V m* 

with q integer. Clearly the the integer level spacing is approximate, as is the degeneracy, 
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FIG. 1: Ground-state density profile for different values of on-site interaction U in a shallow trap 
with u = 0.3. 

because in the N = 5 case discreteness is still manifest. The anticipated dipole and monopole 
(breathing) oscillations are visualized in Fig. 0J where the density evolution in time of a 
superposition of the ground state and each excited state is shown together with the mean 
position and width of the cloud. 

In a trapped BEC, as the coupling grows stronger the frequency of the dipole mode stays 
constant while that of the monopole mode decreases in the Thomas-Fermi limit in one 
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FIG. 2: Quantum fluctuations in the ground-state density, (Sn) = (n 2 ) — (n) 2 , for different values 
of U in the shallow trap, with to = 0.3. 

dimension it is down to E 2 = \/3uj « 1.73a;. In the present system, one cannot hope to see 
this limit, but there is indeed a drop of the second excitation frequency while the lowest one 
initially stays constant. As the coupling strength approaches and exceeds unity, however, 
the discreteness of the system becomes manifest in a ^/-dependence, although weak, of the 
dipole frequency. The dipole and monopole character of the first and second excited state, 
respectively, are unchanged, as is exemplified in Fig. El 
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Lowest lying modes for N=5; co=0.3 



LU 




FIG. 3: Lowest excitation frequencies in units of the trap frequency w in a shallow trap with 
<jj = 0.3. Full lines represent the results for the three lowest excited states obtained by exact 
diagonalization. Dashed line is the outcome of the Bogoliubov approximation. 



In the limit of weak interactions, one expects the Bogoliubov approximation to be valid. 
This approximation can for an inhomogeneous system be effected by linearizing the Gross- 
Pitaevskii equation around its ground-state solution [2]. The latter equation is obtained 
by replacing the field operators in the Hamiltonian, Eq. (JT|) . by their expectation values, 
a r — > z r = (a r ), resulting in the discrete difference equation 



tt\ |2 ^ I \ ^ 2 2 

U\z r \ z r — —{z r+ i + z r -i) + — r z r . 



(7) 



In fact, solving the full Gross-Pitaevskii equation also for the excited states turned out to 
be easier than diagonalizing the linearized Bogoliubov equations, and therefore the former 
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N=5, U/J=0.1 , (o=0.3. First excited state 



Second excited state 





FIG. 4: Time dependence of the density profile for a system with U = 0.1 and u = 0.3 (weak 
trapping and weak interactions). The different sets of bars show the density at successive time 
instances for a system prepared in a superposition of the ground state and an excited state (left 
panel, first excited state; right panel, second excited state) with the amplitudes 0.92 and 0.4, 
respectively. The curves in the lower two panels show the variation of the center-of-mass position 
and mean squared radius with time, showing the dipole and breathing character, respectively, of 
the oscillations. To enhance visibility, the curve for the mean squared radius has been shifted down 
by an amount equal to its time average, defining the plotted quantity as (Sr 2 ) = (^(t)\r 2 \il:(t)} — r 2 
where r 2 is the time average. 

method was chosen. The result thus obtained for the first excited state is included in Fig. 
El and it can be seen that the Bogoliubov approximation fails appreciably even for weak 
interactions. The reason is that the derivation of the Gross-Pitaevskii equation assumes 
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N=5, U/J=100, co=0.3. First excited state 



Second excited state 






-0.2 



E 1 I 



FIG. 5: Same as Fig. 0J but here the coupling is U = 100. 



that the number of atoms in the ground state, N , greatly exceeds the quantum fluctuations 
around it, which, however, are of order unity; for the case of five particles this condition is 
certainly not met. w| 



IV. MODES IN TIGHT TRAPS 

Let us now turn to the case of a tight trap. Figure |H1 displays the ground-state density 
distribution for values of U ranging from weak to strong interactions, in a trap with frequency 
uj = 4.0. The physics is here determined by the balance between interactions and trap and 
the tunneling has little effect. The phase diagram has much more structure in this limit 
compared to the shallow-trap case. 

For weak interactions (U ^$ 0.1), all the particles simply gather in the central well (or 
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FIG. 6: Ground-state density in the case of a tight trap, uj = 4.0. 

the two central wells, if L were even). With increasing U the density distribution flattens, 
and for strong enough interactions the ground state resembles a Mott insulating state, with 
one particle in each of the centermost wells. (If the even/odd parity of the number of wells 
and the number of particles do not match, the edge sites will be partially filled.) Figured 
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FIG. 7: Quantum fluctuations in the ground-state density in a tight trap, to = 4.0. 

verifies that in this state the quantum fluctuations are minimal, which warrants the use of 
the term Mott state. The transition to the Mott state in this five-particle system is in fact 
quite sharp. Also in the limit of small U, the quantum fluctuations become small, but they 
only vanish completely in the limit uo — > oo and the transition is not sharp. 
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Between the two extremes, the shape of the atom cloud is determined by a balance 
between interactions and trap potential. It is quite easy to estimate the crossover values of 
U where the system changes between different types of ground state (we avoid speaking of 
phases for this finite system). If the effective size of the system is R sites, the interaction 
energy Ei and trap energy V scale as 

N 2 U 9 9 

Ei — , V ~ u 2 R 2 . (8) 

R 

Balancing these yields R 3 ~ N 2 U/uj 2 or 

The Mott state, where R = N, thus sets in above U = uj 2 N; inserting the present parameters 
we get U = 80. When the particles gather in the centermost well, R is equal to 1; this 
happens when U = uj 2 /N 2 « 0.5 or smaller. This back-of-the-napkin argument is in almost 
quantitative agreement with the exact numerical findings. 

The lowest few mode frequencies as functions of U when uj = 4.0 are displayed in Fig. |HJ 
(Observe that, in order to emphasize physical interpretation, the frequencies were in Fig. El 
given in units of the trap frequency u, but here it is given in units of the tunneling J). The 
dependence on coupling seen in Figs.EKis seen to have clear consequences also for the mode 
frequencies. In the central- well limit, U ^ 1, the frequencies are easily interpreted since the 
trap determines all the physics, and all the eigenstates of the Hamiltonian are pure Fock 
states or superpositions of degenerate Fock states. The lowest two modes are superpositions 
of the two possible states that result when one particle is removed from the central site 
and put in one of the two neighboring ones, | . . . 00410 . . .) and | . . . 01400 . . .). The limiting 
value of the excitation frequency is equal to u 2 /2 = 8, the excess energy of one particle 
being moved to a neighboring site. The energy of the degenerate third and fourth excited 
states is twice this, since they have two particles in in an off-center well. One may expect 
that among the two lowest excitations, the superposition with a minus sign corresponds to a 
dipole mode and the plus sign corresponds to a breathing mode. Indeed, this is confirmed in 
Fig. IH1 where the time evolutions of the modes are visualized, but note the small amplitude 
of the oscillations. The time dependence is merely due to the small deviations from perfect 
confinement that are still left in the moderately strong trap with frequency u = 4.0; it can 
easily be realized that no current can flow in a system whose eigenstates are pure Fock states 
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U/J 

FIG. 8: Lowest few excitation frequencies in a tight trap, u = 4.0. Dashed lines represent the 
Bogoliubov approximation. 



in configuration space. To see this, consider preparing an initial state as a superposition of 
two eigenstates of the system, — a\A) + (3\B), and study its time evolution, 

\^{t)) = e~ lEAt (a\A) + e-^ ABt (3\B)) , (10) 

where ujab = E B — Ea- Now expand the time evolution in Fock states \f), and obtain for 
the density at the position r, 

(mHa r \m) = E Hf\ A ) + e'^ ABt f3(f\B)\ 2 (f\ata r \f), (11) 
/ 

and we see immediately that if \A) and \B) have no Fock state components in common, 
there can be no time dependence. This fact lies behind the insulating nature of the Mott 
state, but as we have seen, it also prohibits the dynamics in the limit of a very strong trap. 
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N=5, U/J=0.1, co=4. First excited state 



Second excited state 





FIG. 9: Visualization of the two lowest modes for a tight trap, uj = 4.0, and weak coupling, U = 0.1. 
Panels are as in Fig. |1J 

As the coupling gets stronger, the degeneracy is lifted and an intricate pattern of level 
crossings follows; however, the dipole/monopole character of the two lowest modes is retained 
until U approaches the value of approximately 50, where a mixing of the dipole and monopole 
modes starts to be visible, as seen in Fig. EH When U ^ 100, we are in the Mott state and 
the mixing is complete so that the four lowest-lying modes are degenerate with an energy 
equal to 5u 2 /2; this is the energy difference between the second and third well from the 
center as each excited state corresponds just to a displacement of a particle from the center. 
In fact, each of the two lowest excited states has very high overlap (the overlap is 0.98 when 
U = 100) with one of the pure Fock states | . . . 011111000 . . .) and | . . . 000111110 . . .), i. e., 
the ground-state configuration displaced from the center by one lattice site. The third and 
fourth excited states are separated from the first two by an exceedingly small energy gap 
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N=5, U/J=40, co=4. First excited state Second excited state 




1 2 3 1 2 3 

E,t E 2 t 

FIG. 10: Same as Fig. but here the trapping frequency is u = 4.0 and the coupling is U = 40. 

(about 0.005 in units of J), and are mainly superpositions of the state | . . . 010111100 . . .) 
and its mirror reflection. As seen in Fig. ^2 excitation of one of these modes does not result 
in oscillation; this is again a consequence of the fact that the eigenstates are pure Fock states 
in the spatial representation. Indeed, this is a Mott insulating state and excitation does not 
result in particle flow. 

Returning to the mode frequency plot in Fig. |HJ the frequencies have also been calculated 
in the Bogoliubov approximation and are included as dashed lines. Somewhat surprisingly, 
the Bogoliubov approximation performs better in the strong-trapping regime than in the 
weak-trapping regime, although one would naively expect the accuracy to be worse when 
tunneling is suppressed: the quantum fluctuations in local density are small. The appar- 
ent paradox is resolved by noting that in this limit the system is in fact a Bose-Einstein 
condensate: all particles occupy the same state, namely the one confined to the central well. 
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N=5, U/J=100, co=4. First excited state 



Second excited state 
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FIG. 11: Same as Fig. 0J but here the trapping frequency is uj = 4.0 and the coupling is U = 100. 



V. CONCLUSIONS 



The mode frequencies of a trapped boson system in an optical lattice have been studied 
with attention to the dependencies on trap strength and interactions. In the weakly-trapped 
limit, it is shown that the mode frequencies have the usual character of dipole oscillations for 
the lowest-lying mode and breathing oscillations for the next-lowest, and the frequencies are, 
despite the low number of particles and the discreteness, seen to approximately approach 
the result for a harmonically trapped gas when the interactions vanish. Discreteness effects 
set in when the interaction energy is comparable to the tunneling. For sufficiently strong 
trapping, the number fluctuations are quenched, and the dynamics therefore absent, in 
two limits: when the interactions are strong, even the five-particle system displays a quite 
sharp transition to a Mott insulating state. When on the other hand the trapping potential 
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dominates completely, all the particles are trapped in the central well and the amplitude of 
any oscillation goes to zero. It turns out that the Bogoliubov approximation is capable of 
approximating the mode frequencies better for strong trapping than for weak trapping: in 
the weak-trapping case the number of particles is too small for the Bogoliubov approximation 
to work, but a strong trap quenches the fluctuations around the Bose-Einstein condensed 
ground state. 
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